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Self-assembled linear structures like giant cylindrical micelles or discotic molecules 
in solution stacked in flexible columns are systems reminiscent of poly disperse poly- 
mer solutions, ranging from dilute to concentrated solutions as the overaU monomer 
density and/or the chain length increases. These supramolecular polymers have an 
►^ I equilibrium length distribution, the result of a competition between the random 

Q>^ ■ breakage of chains and the fusion of chains to generate longer ones. This scission- 

f^ I recombination mechanism is believed to be responsible of some peculiar dynamical 

Cn i properties like the Maxwell fluid rheological character for entangled micelles. Sim- 

!^ I ulations employing simple mesoscopic models provide a powerful approach to test 

^^ ■ mean-field theories or scaling approaches which have been proposed to rational- 

f^*) I ize the structural and kinetic properties of these soft matter systems. In the present 

work, we review the basic theoretical concepts of these "equilibrium polymers" and 
C^ I some of the important results obtained by simulation approaches. We propose a new 

version of a mesoscopic model in continuous space based on the bead and FENE 
spring polymer model which is treated by Brownian Dynamics and Monte-Carlo 
binding/unbinding reversible changes for adjacent monomers in space, characterized 
Q I by an attempt frequency parameter uj. For a dilute and a moderately semi-dilute 

O i state-points which both correspond to dynamically unentangled regimes, the dy- 

namic properties are found to depend upon u; through the effective life time Tf, of 
the average size chain which, in turn, yields the kinetic reaction coefficients of the 
mean-field kinetic model proposed by Gates. Simple kinetic theories seem to work 
t^ . for times t > Tt, while at shorter time, strong dynamical correlation effects are ob- 

served. Other dynamical properties like the overall monomer diffusion and the mo- 
bility of reactive end-monomers are also investigated. 
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1.2 Introduction 

In the field of self-assembling structures, supramolecular polymers are attracting 
nowadays much attention [1]. It is well known and schematically illustrated in figure 
1 . 1 that some surfactant molecules in solution can self-assemble and form wormlike 
micelles [2]. 

Such micellar solutions exhibit fascinating rheological behaviour, such as shear- 
banding [3], shear-thickenning[4]. Maxwell fluid behaviour [2] or anomalous dif- 
fusion (Levy flight) [5]. Self-assembled stacks of discotic molecules and chains of 
bifunctional molecules are other examples of supramolecular polymers. All these 
examples differ by the nature of the intermolecular forces involved in the self- 
assembling of the basic units, but they lead to a similar physical situation bearing 
much analogy with a traditional system of polydisperse flexibles polymers when 
their length becomes sufficiently large with respect to their persistence length. The 
specificity and originality of these supramolecular polymers comes from the fact that 
these chains are continuously subject to scissions at random places along their con- 
tour and subject to end to end recombinations, leading to a dynamical equilibrium 
between different chain lengths species. These supramolecular polymers are typical 
soft matter systems and the chain length distribution which determines their prop- 
erties, is very sensitive to external conditions (temperature, concentration, external 
fields, salt contents, etc.). 

Different times scales and length scales are involved in these systems, as in 
usual polymer systems. Unlike the latter, which have been much studied during the 
last 50 years[6, 7, 8], theoretical work on micellar solutions are quite recent[l, 2]. 
While building up a molecular scale theory is clearly too heavy for such complex 
systems, theoretical approaches, based on mean-field concept and non-local phe- 
nomenological approaches, can explore with success some equilibrium, and rheolog- 
ical properties[9, 10, 11], under reasonable but often drastic approximations. Con- 
cerning specifically micelle kinetics, recently, O'Shaughnessy and Yu [12] suggested 
that there are two possible kinds of kinetics associated with scission/recombination: 
diffusion controlled and mean-field, which may be distinguished by the time depen- 
dence of the first recombination times distribution function of the chains. 

In parallel, mesoscopic scale computer simulations can shed much light on these 
rich but complex systems, thanks to techniques borrowed from simulations of poly- 
mers. Systems of wormlike micelles can be modelled by "Equilibrium polymers" 
(EP), sometimes called "living polymers", which are polymer chains endowed with 
scission/recombination (S-R) processes taking place in them. The advantage of nu- 
merical studies is their ability to make links either to experiments or to mean- 
field type theories, given that those systems are difficult to characterize by ex- 
periments and the results are not simple to interpret. The bond fluctuation model 
(BFM) has been the object of intensive studies on the statics and dynamics of living 
polymers[13, 14]. 

Brownian dynamics studies of a similar model [15] have been used to study the 
scaling predictions on the dependence of the average chain length upon the overall 
monomer density. In a recent study. Padding and Boek[16] showed that a mesoscopic 
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Fig. 1.1. Some surfactant molecules in solution self-assemble and form long wormlike mi- 
celles which continously break and recombine. Their mass distribution is, hence, in thermal 
equilibrium and they present an important example of the vaste class of systems termed "equi- 
librium polymers"[l].The free energy E of the (spherical) end cap of these micelles has been 
estimated [2] to be of order of lOfcsT. This energy penalty (together with the monomer den- 
sity) determines essentially the static properties and fixes the ratio of the scission and recombi- 
nation rates, ks and kr- Additionally, these rates are influenced by the barrier height B which 
has been estimated to be similar to the end cap energy. Both important energy scales have 
been sketched schematically as a function of a generic reaction coordinate q (see chapter 8 of 
reference [19]). Following closely the analytical description [2, 13] these micellar systems are 
represented in this study by coarse-grained effective potentials in terms of a standard bead- 
spring model. The end cap free energy becomes now an energy penalty for scission events, 
i.e., the creation of two unsaturated chain ends. The dynamical barrier is taken into account by 
means of an attempt frequency uj — exp {—B/kBT). If uj is large, successive breakage and 
recombination events for a given chain can be assumed to be uncorrelated and the recombina- 
tion of a newly created chain ends will be of standard mean-field type. On the other hand, the 
(return) probability that two newly created chain ends recombine immediately must be partic- 
ulary important at large uj. These highly correlated "diffusion controlled" [12] recombination 
events do not contribute do the effective macroscopic reaction rates which determine the dy- 
namics of the system. The key task of this paper is to compare different methods to compute 
the effective rates as a function of the attempt frequency co for a weakly semi-dilute state point 
in the non-entangled dynamical regime. 



model of wormlike micelle known as the FENE-C model[17], seems to obey the 
diffusion-controlled kinetics. 
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In our work[18] which is the basis for this course, we investigate the structural 
and kinetic properties of equilibrium polymers in solution by means of a Brownian 
Dynamics algorithm applied to a standard polymer model[20] coupled to a time- 
reversible EP algorithm for polymer scission-recombinations [13, 14, 21, 22, 23]. 
This particular continuous space model has the advantage that its kinetic properties 
are well separated from the structural and thermodynamical aspects. Our model does 
not allow for cyclic structures and hence, we simplify the computational load and, 
more importantly, the theoretical interpretation of the static and dynamical properties 
of our samples. Anyway, for cylindrical micelles, cyclic structures are expected to 
exist in very low concentrations [2]. 

We study the kinetics of equilibrium polymers solutions at two state points, ex- 
ploring in each case a range of rates of the scission-recombination process. The two 
state points correspond to a dilute solution and a semi-dilute solution regime close 
to the cross-over region. A direct measure of chain overlap can be provided by the 
ratio of the mean distance between chain centres of mass over the radius of giration 
of the average size chain (respectively w 10 and w 50 for the dilute and semi-dilute 
cases). This gives ratios of 1.4 for the dilute case and 0.6 for the semi-dilute case, 
which suggests that even in the semi-dilute case, the chain relaxation is still typically 
Rouse like and thus yet kinetically entangled. 

The main focus is to link the microscopic model to the more macroscopic kinetics 
theories, e.g. by Gates [2], and to discuss in a well defined case the concepts of 
'diffusion controlled' and 'mean-field', defined by authors of[12, 16]. In particular, 
we want to estimate the effective reaction rates for the scission and recombination 
processes comparing different methods. Furthermore, we shall study the diffusion of 
monomers and discuss their link to the scission/recombination kinetics. 



1.3 Theoretical framework 

1.3.1 Statistical mechanics derivation of the distribution of chain lengths 

To treat a system of supramolecular polymers theoretically, it is convenient to work at 
the mesoscopic scale using a model of linear flexible polymers made of L monomers 
of size b linked together by a non permanent bonding scheme. Within the system, 
individual chain lengths fluctuate by bond scission and by fusion of two chain ends 
of different chains. Statistical mechanics can be employed to predict the equilibrium 
distribution of chain lengths[2]. In terms of the equilibrium chain number density 
co{L), the average chain length Lq and the total monomer density cj) are given by 

^0 — v^oo 7TT- K^-^) 

°° M 

0=^ico(L) = -- (1.2) 

L=0 

where M denotes the total number of monomers in the system. 
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Conceptually, we consider the Helmholtz free energy F{V, T; {N{L)}, Ng) of a 
mixture of chains molecules of different length L in solvent where, in addition to the 
temperature T, the volume V and the number of solvent molecules Ng, the number of 
chains of each specific length N{L) is fixed. Let F{V, T; M, Ng) be the Helmholtz 
free energy of a similar system where only the total number of solute monomers M 
is fixed. The equilibrium chain length distribution co{L) = N{L)/V will result from 
the set {N{L)} which satisfies the condition 



FiV,T;M,Ng)=rmn^^NiL)} 



FiV, T; {N{L)}, TV.) + M ^Z LN{L) 

L=0 



(1.3) 



The parameter /i is the Lagrange multiplier associated with the constraint that 
individual numbers of chains N{L) must keep fixed the total number of monomers 
M = J2T LN{L). Minimisation requires that the first derivative with respect to any 
N{L') variable (L'=l,2, ...) is zero, giving 

We expect the entropic part of the total free energy F{V, T; {N{L)}, Ng) to be 
the sum of translational and chain internal configurational contributions which both 
depend upon the way the M monomers are arranged into a particular chain size dis- 
tribution. For the translation part, the polydisperse system entropy is estimated as the 
ideal mixture entropy Sid 

SUV, T; {N{L)}, Ng) ^ -fcs ^ N{L) In {CN{L)) + 5^°'" (1.5) 

L 

where C — b^/V is a dimensionless constant independent of L and where 5''*°'" is 
the solvent contribution, independent of the N{L) distribution. The configurational 
entropy of an individual chain with L monomers is written as Si{L) = fcs In Hl, in 
terms of Hl, the total number of configurations of the chain. Adding the configura- 
tional contributions to Sid as given by eq.(1.5), the total entropy becomes 

SiV, T; {N{L)}, Ng) = -kg J2 ^W P^ iCN{L)) - In r^i] + S'°'" (1.6) 

L 

We now turn to the energy E{V, T; {N{L)}, Ng) of the same system. If Ei{V, T; 
L) represents the internal energy of a chain of L monomers and Eg the energy of a 
solvent molecule, the energy can be written as 

^(F,r;{iV(i)},iV,)=^A^(L)^i(F,T;L)+iV,^,(F,T) (L7) 



The key contribution in Ei is the chain end-cap energy E which corresponds to 
the chain end energy penalty required to break a chain in two pieces. We will suppose 
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that Ei{L) ~ E + Le where e is an irrelevant energy per monomer as Me, its total 
contribution to the system energy, is independent of the chain length distribution. 

The present approximation of the total free energy of the system is thus given by 
incorporating in the general expression (1.3) the expressions (1.6) and (1.7), giving 

PF(y, T; {N{L)}, Ns) = ^ N{L) [In N{L) + In C - In I2l + I^E + (31L] 

(1.8) 
where irrelevant constant solvent terms have been omitted as we only need the first 
derivative of the free energy with respect to N{L), which now takes the form 

^ =lnN(L) + \nC -InilL+pE + piL + l. (1.9) 



SN{L) 
With this expression, the minimisation condition on N(L) becomes 

IniV(i) + InC - In Ql + PE + 1 + fi' L = (1.10) 

where /i' = {/3n + (3e). We note at this stage that the second derivative of l3F{V, T; 
{N{L)}) + /x' J2l LN{L) with respect to N{L) and N{L') variables gives the non 
negative result j^tt, indicating that the extremum is indeed a minimum. 
Solving for N{L) in eq.(l.lO), we get 

N{L) ^ C'^^ exp -ifi' L + f3E - \n Hl) (1.11) 

where C" — eC while, according to eq. (1.2)), /i' must be such that 

^Lexp-(/i'i + /3£-lnl2L) = MC' (1.12) 

L 

The equilibrium N{L) variables are also related to the equilibrium chain length 

average Lq (see eq. (1.1)), so that 

MC 
yexp-(^'L + /3£;-lnr2L)==— — (1.13) 

L ^0 

To progress, we now need to specify the explicit L dependence of il^. The tradi- 
tional single chain theories of polymer physics provide universal expression of Hl in 
terms of the polymer size, the environment being simply taken into account through 
the solvent quality and the swollen blob size in the semi-dilute (good solvent) case. 

The case of mean field or ideal chains 

The basic mean-field or ideal chain model for a L segments chain gives 

nf^lCiz'^] (1.14) 
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where z is the single monomer partition function and C\ a dimensionless constant. 
Adapting eq.(l.ll), one has 

iV(L) = ^exp-(/3S)cxp(-A^"L) (1.15) 

where /i" = /i' — In z must, according to eq. (1.12), be such that 

1 MC 

while eq.(1.13) takes the form 

^exp-(M"L) = \= , ^ ^'^' ,.^, (1.17) 

In eqs.(1.16) and (1.17), sums over L from 1 to oo have been approximated by 
the result of their continuous integral counterparts. 

Combining eqs.(1.15), (1.16) and (1.17), one gets the final expression for the 
chain number densities 

co(i) = 4exp(-^) (1.18) 



^0 

with the average polymer length given by 



(3E 



Lo = Bi/20^cxp(^^j. (1.19) 

where B = eb'^ /Ci is a constant depending upon the monomer size b and the pref- 
actor in the number of ideal chain configurations in eq.(1.14). 

The case of dilute chains in good solvent 

Polymer solutions in good solvent are in a dilute regime when chains do not overlap 
and in semi-dilute regime when chains do strongly overlap while the total monomer 
volume fraction is still well below its melt value. In the semi-dilute regime, chains 
remain swollen locally over some correlation length, known as the swollen blob size 
X, but they are ideal over larger distances as a result of the screening of excluded 
volume interactions between blobs. 

Specifically, for a given monomer number density (f>, the blob size is given by the 
condition that the blob volume times must be equal to the number of monomers 
L* in the swollen blob. This gives in terms of the reduced number density (/>' — b^cf) 

L* - (/,'(t^) (1.20) 

X-60'(t^) (1.21) 

where u = 0.588 in present good solvent conditions [8]. 
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In living polymers characterized by a monomer number density (j> and some av- 
eraged chain length L^, the semi-dilute conditons correspond to the case Lq >> L*. 
We discuss in this subsection the theory for the dilute case where L^ << L*. We 
will come back to the semi-dilute case in the next subsection. 

Self-avoiding walks statistics apply to dilute chains in good solvent, and we thus 
adopt the number of configurations [6, 7](See especially page 128 of the book of 
Grosberg and Khokhlov[7]) 

r2f ^ = CiL^'^-^h^ (1.22) 

where 7 is the (entropy related) universal exponent equal to 1.165. 
Incorporating expression (1.22) in eq. (1.11), one gets 

N{L) = -^cxp-(/3£;)i('^^i)exp-(^"L) (1.23) 

where B was introduced in eq.(1.19) and where /i" = /_/ — Inz must be fixed by 
eq.(1.2) 

^L^exp-(/^"i) =B0exp(/3^) (1.24) 

L 

while eq.(1.13) takes here the form 

y L('^-i'exp-(^"i) = ^exp(/3£;) (1.25) 

L ^" 

If L is treated as a continuous variable, eqs(1.24) and (1.25) can be rewritten in 
terms of the Euler Gamma function satisfying r{x) = xr{x — 1) as 

/ L''cxp^{n"L)dL = „|^^^^^ = S0cxp (/?£:) (1.26) 

/"L(^-i)cxp-(M"L)di = ^ - ^exp(/3i?) (1.27) 

Jo M ' -f^o 

From eqs (1.26) and (1.27), one gets 

m" = ^ (1-28) 

i30cxp(/3i?) = :^;il±il4''+i) (1.29) 

These results lead then finally to the Schulz-Zimm distribution of chain lengths, 
namely 

coiL) = '-^^^^Li^~^^e.pi-,^) (1.30) 

and an average polymer length given by 
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The semi-dilute case 

We consider here the semi-dilute case in good solvent where the average length of 
living polymers Lq is much larger than the blob length L*. The usual picture of a 
semi-dilute polymer solution is an assembly of ideal chains made of blobs of size x- 
Using this approach. Gates and Candau [2] and later J. P. Wittmer et al [13] derived 
the relevant equilibrium polymer size distribution. In this subsection, we adapt their 
derivation to the theoretical framework presented above. 

Let ill, be the number of internal configurations per blob and z' some coordina- 
tion number for successive blobs. As there are n^ = L/L* blobs for a chain of L 
monomers, we write the total number of internal configurations of a chain of size L 
as 

nl" ^ Cii*(T-i)|2,''/^*z'^/^* (1.32) 

where 7 is the universal exponent in the excluded volume chain statistics met earlier 
for chains in dilute solutions. The important factor i*(''^i) can be seen as an entropy 
correction for chain ends just like E was an energy correction to Li. This entropic 
term which involves the number of monomers per blob, is needed to take into account 
that when a chain breaks, its two ends are subject to a reduced excluded volume 
repulsion. The other factors in eq.(1.32) will lead to terms linear in L after taking the 
logarithm and thus will be absorbed in the Lagrange multiplier definition, as seen 
earlier in similar cases for ideal and dilute chains. The resulting expression of N{L) 
in terms of the Lagrange multiplier (cfr eq.(1.15)) can then be written by analogy as 

N{L) = ^ exp - [I3E - (7 - 1) hiL*] exp (-/iL) (1.33) 

Proceeding as in the ideal case (simply replacing at every step the constant (3E 
by f3E — (7 — 1) In L*, one recovers in the semi-dilute case the simple exponential 
distribution 

co(L) = -^exp(-^) (1.34) 

with a slightly different formula for the average polymer length 

Lo ex r exp ('^'j (1.35) 

where a = ^(1 + 3^" \ ) is about 0.6. 

1.3.2 A kinetic model for scissions and recombinations 

The interest for wormlike micelles dynamics came from the experimental observa- 
tion that entangled flexible supramolecular polymers display, after an initial strain, 
a simple exponential stress relaxation which is qualitatively different from the be- 
haviour of usual entangled melts. In the latter system, individual chains must leave 
by a reptation mechanism the strained topological tube created by the entangled 
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temporarily network, in order to relax the shearing forces. A plausible scission- 
recombination model, allowing for such an extra relaxation mechanism, was shown 
[2] to lead to an exponential decay of the shearing forces with a decay Maxwell time 
equal to r w jT\f^ where r^ep is the chain reptation time and t\^ is the mean life 
time of a chain of average size in the system. 

We will assume in the following that the Gate's scission-recombination model 
governing the population dynamics, originally deviced to explain entangled equilib- 
rium polymer melt rheology, should also apply to the kinetically unentangled regime 
which is explored in the present work. 

This kinetic model [9] assumes that 

• the scission of a chain is a unimolecular process, which occurs with equal prob- 
ability per unit time and per unit length on all chains. The rate of this reaction is 
a constant kg for each chemical bond, giving 

n = 7^- (1-36) 

for the lifetime of a chain of mean length Lq before it breaks into two pieces. 

• recombination is a bimolecular process, with a rate fc^ which is identical for all 
chain ends, independent of the molecular weight of the two reacting species they 
belong to. It is assumed that recombination takes place with a new partner respect 
to its previous dissociation as chain end spatial correlations are neglected within 
the present mean field theory approach. It results from detailed balance that the 
mean life time of a chain end is also equal to T(,. 

Let c(t, V) be the number of chains per unit volume having a size L at time t. 
Gates [9] writes the kinetic equations as 



^^^hllL ^ -k,Lc(t, L) + 2k, I c(t, L')dL' (1.37) 

dt Jl 

k f^ f°° 

+ -^ c{t, L')c{t, L - L')dL' - krc{t, L) / c(t, L')dL'(\3%) 
2 Jo Jo 

where the two first terms deal with chain scission (respectively disappearence or 
appearence of chains with length L) while the two latter terms deal with chain re- 
combination (respectively provoking the appearence or disappearence of chains of 
length L). 

It is remarkable that the solution of this empirical kinetic model leads to an ex- 
ponential distribution of chain lengths. Indeed, direct substitution of solution cq{L) 
in the above equation leads to the detailed balance condition: 

ct>^^=Ll (1.39) 

the ratio of the two kinetic constant being thus restricted by the thermodynamic state. 
Detailed balance means that for the equilibrium distribution cq{L), the number of 
scissions is equal to the number of recombinations. The total number of scissions 
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and recombinations per unit volume and per unit time, denoted respectively as n^ 
and Tir, can be expressed as 



Us 


— '^s j2 
~ 9 r4 


1 L c-np {—L /LQ)dL ~ ks<p 
'o 

/ dL' / <ii"exp(- — )exp(- 
Jq Jq Lq 


-£'^ 


02 fc. 


(1.40) 
(1.41) 


Tlr 


- 2Ll 



and it can be easily verified that detailed balance condition implies n^ = rir. 

Mean field theory assumes that a polymer of length L will break on average 
after a time equal to Tb = {ksL)~-^ through a Poisson process. This implies that the 
distribution of first breaking times (equal to the survival times distribution) must be 
of the form 

if'(t) = exp(- — ) (1.42) 

n 

for a chain of average size. Detailed balance then requires that the same distribution 
represents the distribution of first recombination for a chain end[2]. Accordingly, 
throughout the rest of this chapter, the symbol Tb will represent as well the aver- 
age time to break a polymer of average size or the average time between end chain 
recombinations. In the same spirit, we stress that among the different estimates of 
Tb proposed in this work, some are based on analyzing the scission statistics while 
others are based on the recombination statistics. 

Two additional points may be stressed at this stage: 

• The mean field model in the present context has been questionned [12] because 
in many applications, there are indications that a newly created chain end often 
recombines after a short diffusive walk with its original partner. In that case, a 
possibly large number of breaking events are just not effective and the kinetics 
proceeds thus differently. 

• Given the statistical mechanics analysis in the previous subsection, we see that 
the equilibrium distribution of chain lengths resulting from the simple empirical 
kinetic model is perfectly compatible with the equilibrium distribution in poly- 
mer solutions at the 9 point (ideal chains) or for semi-solutions (ideal chains of 
swollen blobs). The kinetic model is still pertinent, given its simplicity, in the 
case of dilute solutions as the static chain length distribution based on statistical 
mechanics, although non exponential, is not very far from it. 

1.4 The mesoscopic model and its Brownian Dynamics 
implementation 

In order to test the theoretical understanding summarized in the previous section, 
we have devised a new microscopic model bearing resemblances with an off-lattice 
Monte-Carlo model proposed earlier [21], but which uses Brownian Dynamics to 
follow the monomer dynamics for a standard polymer chain model[20]. To be able 
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to clearly unravel static and dynamic aspects, our model contains a kinetic param- 
eter which allows to vary the rates of the reactions without affecting the static and 
thermodynamic properties. As already mentioned in the introduction, two thermody- 
namic states corresponding to equilibrium polymers dissolved in a good solvent are 
considered, namely a dilute and a semi-dilute solution. 

We consider a set of micelles consisting of (non-cyclic) linear assemblies of 
Brownian particules. Within such a linear assembly, the bounding potential Ui{t) 
acting between adjacent particles is expressed as the sum of a repulsive Lennard- 
Jones (shifted and truncated at its minimum) and an attractive part of the FENE 
type[20]. The pair potential U2{r) governing the interactions between any un- 
bounded pair (both intramicellar and intermicellar) is a pure repulsive potential cor- 
responding to a simple Lennard- Jones potential shifted and truncated at its minimum. 
This choice of effective interaction between monomers implies good solvent condi- 
tions. 
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Fig. 1.2. Bounded potential Ui (r) (continuous curve) and unbounded potential f/2 (r) (dashed 
line) between a pair of monomers. Each monomer can participate to at most two bounded 
links (with potential Ui{r)). All unbounded monomer pairs interact via the potential U2{r). 
E is a parameter tuning the energy required to open the bond. The figure also shows the F re- 
gion where potential swaps (equivalent to bond scissions or bond recombinations) are allowed 
during the brownian dynamics simulation is indicated. 
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Using the Heaviside function 0{x) = or 1 for a; < or a; > respectively, 
explicit expressions (see figure 1.2) are 



r r 4 



0{2^/^a-r)) (1.43) 



U2ir) = 4e 

Uiir) = U2{r) - Q.hkR' In [l - [-]^ j - t/,„,„ - E (1.44) 

In the second expression valid for r < R, k — 30e/cr^ is the spring constant and 
R = 1.50- is the value at which the FENE potential diverges. Umin is the minimum 
value of the sum of the two first terms of the second expression (occuring rmin = 
0.96094(7 for the adopted parameters) while E isa key parameter which corresponds, 
when pair potential exchanges between bounded and unbounded situations will be 
allowed, to the typical energy gain (loss) when an unbounded (bounded) pair is the 
object of a recombination (scission). 

We deal with M monomers at temperature T enclosed in a volume V, defining a 
monomer number concentration cj) = M/V. The relevant canonical partition func- 
tion is a sum over all microscopic states specified not only by monomer positions 
(and momenta) but also by a set of variables defining unambiguously the connectivity 
scheme for the particular microscopic state. Disregarding the irrelevant momentum 
variables in Brownian Dynamics, we write the partition function as 

Qmvt-EE°-p(-^^t¥^) (1.45) 

{r} m ^ 

where the double sum runs over all possible arrangements of M distinguishable 
monomers and over all related bounding schemes. The first sum runs over all spatial 
positions {r} of the M monomers, each configuration being characterized by a list 
of K monomer-monomer pairs with relative distance shorter than R. The second sum 
runs, for a given spatial arrangement with a subset of K pairs available for bound- 
ing, over all allowed distinct bounding schemes, namely those among the 2^ distinct 
possibilities which satisfy two restricting rules, namely 

• no monomer can be engaged into more than two bounding pairs (= no branching). 

• no cyclic bounding structure is allowed. 

The potential energy associated to a particular configuration can be written 

U{{r}.m) = l^^(s,,[/i(r,,) + (1 - s,,)U2{n,)) (1.46) 

where Sij = for all pairs with Vij > R while Sij = either or 1 for the K pairs 
with rij < R. Index m in the two previous expressions corresponds to a particular 
bounding network representing a particular set of and 1 for the K Sij values. 

In order to discuss topological changes (change in bounding network) during the 
phase space exploration, it is useful to consider that each monomer possesses two 
arms available for bounding. If both arms are free, we have an isolated monomer 
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with degree of polymerisation L ~ 1. If at least one arm of a particular monomer is 
linked through potential Ui{r) to an arm of another monomer, these two monomers 
belong to the same linear structure (micelle) of degree of polymerisation L > 2. A 
terminal or an interior monomer will correspond to a monomer engaged respectively 
in one or two bonds with other monomers. 

The phase space exploration will be performed by a Brownian Dynamics (BD) 
scheme (controling r positions updates) coupled to a single bond creation (replace- 
ment of potential U2{r) by Ui(r)) I single bond annihilation (replacement of poten- 
tial Ui{r) by U2{r)) scheme satifying microscopic reversibility. 

The BD is standard and samples the set of positions of the M monomers accord- 
ing to a canonical ensemble corresponding to a particular potential energy (equiva- 
lent in the present case to a particular bonding network m). 

The single BD step of particule i subject to a total force Fi is simply 

n{t + At) = nit) + ^At + R,{At) (1 .47) 

where the last term corresponds to a vectorial random gaussian quantity with first 
and second moments given by 

<RiaiAt)> = (1.48) 

< RU^t)R,p{At) > = 2(-^)Z\i,5„^5y (1.49) 

for arbitrary particles i and j and where af3 stand for the x,y,z cartesian components. 
At this stage, it is useful to fix units. In the following, we will adopt the size of 
the monomer a as unit of length, the e parameter as energy unit and we will adopt 
^a'^ /(Sne) as time unit. We also introduce the reduced temperature ksT/e = T*, so 
that in reduced coordinates (written with symbol*) the algorithm becomes 

<R*^{At*)>^0 (1.50) 

< R*JAt*)R*p{At*) > = l-At*Sa,p6,j T* (1.51) 

which means that each particle, if isolated in the solvent, would diffuse with a RMSD 
of ^2T* j-K per unit of time. 

The bonding network is itself the object of random instantaneous changes pro- 
vided by a Monte-Carlo algorithm which is built according to the standard Metropo- 
lis scheme. 

The probability Pra,n to go from a bounding network m to a different one n is 
written as 

p \trial Tjacc /i co-i 



where A^*°' is the trial probability to reach a new bounding network n starting 
from the old one m, within a single MC step. This trial probability is chosen here 
to be symmetric as usually adopted in Metropolis Monte Carlo schemes. This trial 
probability is chosen to be different from zero only if both bounding networks n and 
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m differ by the status of a single bond, say the pair of monomers (ij). To satisfy 
microreversibilty, the acceptance probability P^^ for the trial (rn — > n) must be 
given by 

^=^M»[l.oxp(-MlLjtMiH),|. ,,33^ 

Kb-L 

In the present case, as a single pair (ij) changes its status, the acceptance proba- 
bility takes the explicit form 

PZn = M.n[l,exp(.(M:4-i^lM))] (1.54) 

for bond scission and bond recombination respectively. 

The way to specify the trial matrix A*^*°' starts by defining a range of distances, 
called r and defined by 0.96 < r < 1.20 within the range r < R. For r G P, a 
change of bounding is allowed as long as the two restricting rules stated above are 
respected. Consider the particular configuration illustrated by figure 1.3 where M=7 
monomers located at the shown positions, are characterized by a connecting scheme 
made explicit by representing a bounding potential by a continuous line. The dashed 
line between monomers 5 and 6 represents the changing pair with distance r G P 
which is a bond Ui (r) in configuration m but is just an ordinary intermolecular pair 
U2{r) in configuration n. 

For further purposes, we have also indicated in figure 1.3 by a dotted line all 
pairs with r E P which are potentially able to undergo a change from a non bounded 
state to a bounded one in the case where the (56) pair, on which we focus, is non 
bounded (state n). Note that bond (35) is not represented by a dotted line eventhough 
the distance is within the P range: a bond formation in that case is not allowed as it 
would lead to a cyclic conformation. Note also that in state n, monomer 5 could thus 
form bonds either with monomer 2 or with monomer 6 while monomer 6 can only 
form a bond with monomer 5. 

We now state the algorithm and come back later on the special (m — > n) tran- 
sition illustrated in figure 1.3. 

During the BD dynamics, with an attempt frequency ui per arm and per unit of 
time, a change of the chosen arm status (bounded Ui{r) to unbounded U2{r) or 
unbounded to bounded) is tried. If it is accepted as a "trial move" of the bounding 
network, it obviously implies the modification of the status of a paired arm belonging 
to another monomer situated at a distance r € P from the monomer chosen in the 
first place. 

The trial move goes as follows: a particular arm is chosen, say arm 1 of monomer 
i, and one first checks whether this arm is engaged in a bounding pair or not. 

• If the selected arm is bounded to another arm (say arm 2 of monomer j) and 
the distance between the two monomers lies within the interval r,;j e P, an 
opening is attempted with a probability l/{Ni + 1), where the integer Ni rep- 
resents the number of monomers available for bonding with monomer i, besides 
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the monomer j {Ni is thus the number of monomers with at least one free arm 
whose distance to the monomer carrying the originally selected arm i lies within 
the interval F, excluding from counting the monomer j and any particular arm 
leading to a ring closure). If the trial change consisting in opening the (ij) pair 
is refused (either because the distance is not within the F range or because the 
opening attempt has failed in the case Ni > 0), the MC step is stopped without 
bonding network change (This implies that the BD restarts with the (ij) pair being 
bonded as before). 
• If the selected arm (again arm 1 of monomer i) is free from bonding, a search is 
made to detect all monomers with at least one arm free which lie in the "reactive" 
distance range r £ F from the selected monomer i (Note that if monomer i is a 
terminal monomer of a chain, one needs to eliminate from the list if needed, 
the other terminal monomer of the same chain in order to avoid cyclic mi- 
celles configurations). Among the monomers of this "reactive" neighbour list, 
one monomer is then selected at random with equal probability to provide an 
explicit trial bonding attempt between monomer i and the particular monomer 
chosen from the list. Note that if the Ust is empty, it means that the trial attempt 
to create a new bond involving arm 1 of monomer i has failed and no change in 
the bonding network will take place. 

In both cases, if a trial change is proposed, it will be accepted with the probability 
P^ji defined earlier If the change is accepted, BD will be pursued with the new 
bonding scheme (state n) while if the trial move is finally rejected, BD restarts with 
the original bonding scheme corresponfing to state m. 

Coming back to figure 1 .3, we now show that the MC algorithm mentioned above 
garantees that the matrix A^"' is symmetric, an important issue as it leads to the 
micro-reversibility property when combined with the acceptance probabilities de- 
scribed earlier. Let us define as Farm ~ 1/2N the probability to select a particular 
arm, a uniform quantity. 

If configuration m with pair (56) being "bounded" is taken as the starting config- 
uration, the number of available arms to form alternative bonds with monomer 5 and 
monomer 6 are respectively N^ — 1 and Nq = 0. Therefore, applying the MC rules 
described above, the probability to get configuration n where the pair (56) has to 
be unbounded is given by the sum of probabilities to arrive at this situation through 
selection of the arm of monomer 5 engaged in the bond with monomer 6 or through 
selection of the arm of monomer 6 engaged in a bond with monomer 5. This gives 

1 13 

^m,n = Parm * ^^ . -, + Parm ^ , ^ 9 * "'"'" (1.56) 

If configuration n with pair (56) being "unbounded" is taken as the starting con- 
figuration, the application of the MC rules lead to the probability to get configuration 
m where the pair (56) has to be bounded is given by the sum of probabilities to ar- 
rive at this situation through selection of the free arm of monomer 5 (which has two 
bounding possibilities, namely with monomers 2 and 6) or through selection of the 
free arm of monomer 6 which can only form a bond with monomer 5. This gives 
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Fig. 1.3. Exemplary configuration of a 7 monomers system in state n where monomers 3,4 and 
5 form a trimer and monomers 6 and 7 a dimer. All pairs of monomers with mutual distances 
within the F region are indicated by a dotted or a dashed line. In the text, we consider the 
Monte Carlo scheme for transitions between states n and m which only differ by the fact 
that in state n and m the 5-6 pair is respectively open or bounded. The n — » m transition 
corresponds to the creation of a pentamer by connecting a dimer and a trimer while the m — > 
n transition leads to the opposite scission. The cross symbol on link 3-5 indicates that in state 
n, when looking to all monomers which could form a new link with monomer 5, monomer 3 
is excluded because it would lead to a cyclic polymer which is not allowed within the present 
model. 



A. 



trial 



1 3 



(1.57) 



showing the required matrix symmetry. 



1.5 Characterization of static properties 

In the following, all quantities will systematically be expressed in reduced units with- 
out an explicit * symbol. 

1.5.1 List of experiments and cliain length distributions 



We consider a system of M = 1000 monomers interacting via the potential model 
defined by eqs.(1.44) at T = 1. This is a rather small system by nowadays standards 
but we were only interested here to a preliminary demonstrative study before going 
in the future to larger systems and a more efficient description of solvent effects. 
This system is studied at two state points: 

• A dilute solution at the number density = 0.05 and an energy parameter iJ = 8. 

• A semi-dilute solution at the number density (^ = 0.15 and an energy parameter 

E = 10. 
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Each system evolves according to the Brownian Dynamics algorithm with time step 
At = 0.001 and is subject to random trials of bond scission/recombination with arm 
trial attempts frequencies in the range 0.1 < w < 5. This choice was motivated by 
the necessity to speed up the kinetics and hence, getting reasonably good statistics 
on dynamical properties. Obviously, all static properties are independent of ui and 
therefore, all data can be cumulated over all cu values. 

The systems have been equilibrated for a time of 2 10^. Production runs have 
been followed for a similar time Tmn = 2 lO"'. The equilibrium values of the 
mean chain length io and the total number of chains < Nch > are related by 
M = Lo< Nch >■ 

The dilute case 

The first experiment deals with a dilute solution as the average chain length Lq = 
10.4 ± 0.1 is much smaller than the crossover value at that monomer number density 
as calculated by eq.(1.21), L* = 50.5. 



'I I' 



'I I' 



dilute 
semi-dilute 



^^ ^ ^H.'^^^^ j...^'S£cy^ 




Fig. 1.4. Chain end pair distribution function gee{r) in the dilute and semi-dilute experiments. 
Note that the F region where the bounding changes take place correspond to the region of first 
(fast) increase of the distribution function around r = a. 
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Fig. 1.5. Normalized distribution function P{r) of the distance between bounded monomers 
in the semi-dilute experiment. The corresponding function in the dilute solution experiment is 
marginally different from the semi-dilute case and is therefore not shown in the figure. 



The populations of bonds ready to open within the F range is < Nir >= 551.6 
and represents 61% of the bonded pairs < TVi >, while the population of free arms 
pairs ready to close, again within the same F range, is only < N2r >— 1.16. Fig- 
ures 1.4 and 1.5 show respectively the pair correlation function for chain ends (un- 
saturated arms) and the distribution of the distance between bounded monomers, in 
particular within the F region where potential changes do occur. The "free arm" 
fraction, namely 1 - Ni/M = 0.0957, is close to L„^ =< Nch/M >= 0.0964. 

Dilute solution conditions are confirmed by the chain length distribution shown 
in figure 1.6. For dilute conditions, a distribution given by eq.(1.30) is expected. 
Accordingly, we have fitted our data with a single parameter (prefactor) fit function, 
c{L) = AoL'^~^ exp (— 7L/L0) where Lo was given its computed average value and 
where 7 was given its expected value, 7 = 1.165. This curve is significantly better 
than the simple exponential distribution expected for ideal or semi dilute chains. If 
7 is left as a second free parameter in the fit, it takes an even larger value 1.187 and 
the fit (not shown) is shghtly better especially at high L. 
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L/L„ 



Fig. 1.6. Distribution of chain length number densities for the dilute case in good solvent. The 
data (squares) are fitted using eq.(1.30) with imposed value 7 — 1.165. The dashed line shows 
the exponential distribution Co{L) oc exp ( — ■^) which does not fit data as well. 



The semi-dilute case 



In the second state point experiments, the average chain length is found to be 
< L >~ 49.9 ± 0.1, a value which is much larger than the crossover value L* « 12 
at = 0.15 according to eq.(1.20). 

The populations of bonds ready to open in the F region is < Nip >= 597 ± 
1., while the population of free arms pairs ready to close is < N2r >= 0.183. 
Distance distributions of unbounded and bounded pairs are given in figures 1 .4 and 
1 .5 together with the dilute case data. We note that 5(7-) function between free ends is 
not very different between the two cases. The "free arm" fraction l — Ni/M = 0.021 
is 5 times smaller than in the dilute case as the chains are 5 times longer. 

Semi-dilute solution conditions are confirmed by the observation of an exponen- 
tial distribution of the reduced chain length (see figure 1.7). In the latter figure, we 
show the best fit by a one parameter fitting function c^"^* [L) ~ A exp {—L/Lq). 
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L/L„ 



Fig. 1.7. Distribution of chain lengths number densities for the semi-dilute case. The data 
(circles) are fitted with a simple exponential function Aexp{—L/Lo) (continuous curve), 
in reference to eq.(1.34) with imposed average chain size Lo — 49.9. Dashed lines show a 
function exp (- 



-1.165t^) which does not fit the data as well. 



1.5.2 Chain length conformational analysis 



Figures 1.8 and 1.9 show for dilute and semi-dilute state points respectively, the 
average of the square of the radius of gyration and the average of the square of the 
end-to-end vector as a function of the chain length within the polydisperse sample. 

In figure 1.8 relative to the dilute case, standard power law scaling 1?^ law with 
V = 0.588 is indicated for long chains. For our polydisperse system with an average 
chain length of Lq w 10, our data on chains with L > Lq aie affected by too large 
statistical error bars to test critically the validity of this asymptotic regime and its 
exact range of application, limited in the low L regime by usual short chain effects 
clearly visible below L ~ 25. 

The conformational properties of chains in the semi-dilute system are compatible 
with ideal chains for L > L* « 12 where L* is estimated from eq.(1.20) with 
(p = 0.15, and with a larger exponent power law for lower L values but poor statistics 
and small chain effects prevent a firm test of scaling laws for these semi-dilute chains. 

When < R^ > data from both state points are plotted versus L on a unique graph 
(not shown), the short chains behaviours are found to be undistinguishable, while di- 
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Fig. 1.8. Conformational properties in the dilute case: < R > (circles) and < Rg > (squares) 
versus chain length L. The fitting functions which are power law with exponent fixed to 2^ = 
1.176 (prefactors 0.876 for < R^ > and 5.854 for < R^ >) fit the data for longer chains 
only 



vergences between the two state points are observed for long chains. The power laws 
shown in figures 1.8 and 1 .9 intersect at the blob length value of L* = 37 monomers, 
giving a blob size or a correlation length of x =< R^ >^/^~ 20. Note that in this 
analysis, we implictly assume that the prefactor of the dilute chains power law is 
independent of (j> while we consider the (p dependent prefactor in the ideal chain 
power law corresponding to the semi-dilute state at = 0.15. We thus see that the 
blob size as determined from the intersection of asymptotic power laws (neglecting 
short chains effects) suggest a numerical prefactor of about 3 in the scaling relations 
(1.20). On this basis, our dilute system at ^ = 0.05 would be characterized by a 
blob length of L* = 50 and thus a ratio L^/L* = 0.2 suggesting a rather clear 
dilute solution character. The same ratio is Lq/L* = 1.4 for the more concentrated 
solution showing that the usual concept of chains of blobs in a semi-dilute solution 
applies only marginally to the longer chains in the sample. To conclude, this analysis 
indicates that, even in the semi-dilute case, the chain dynamics should not be sig- 
nificantly affected by the weak entanglements between chains. Let us stress that the 
semi-dilute character of the sample is however well illustrated by static properties 
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Fig. 1.9. Conformational properties in the semi-dilute case: < R^ > (squares) and < R^ > 
(circles) versus chain length L. The fitting functions (linear in L) assume ideal chain statistics 
with prefactors 1.862 for < 7?g > and 11.06 for < R^ >. The fit works for the L > L* 
behaviour. 



such as the exponential chain length distribution and the ideal chain power laws for 
chain sizes. 



1.6 Kinetics analysis 



The number of Monte-Carlo "accepted" scissions (or "accepted" recombinations) are 
obtained by simple counting during the simulations. In table 1.1 (dilute case) and in 
table 1.2 (semi-dilute case), we list the number Ug of "accepted" transitions per unit 
time and unit volume, for the two state points and for the different attempt frequen- 
cies Lu investigated. For each state point, we observe that the number of scissions 
per unit of time is, as expected, proportional to w and we have also verified that the 
number of recombinations differs from the number of scissions by marginal amounts 
(0.05%), which shows that the chain length distributions are well equilibrated. 

The quantity a mentioned in the tables is the fraction of scissions which lead 
to a recombination with a new partner. In our analysis, the fate of a scission has 
been initially followed for a maximum time of T^ax ~ 1000. The value quoted in 
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Table 1.1. Kinetic data for the dilute case at four scission-recombination attempt frequencies 
uj. Us is the total number of Monte-Carlo scissions (recombinations) per unit time and per 
unit volume, a is the fraction of scissions where newly created free arms recombine with a 
new partner. Superscript (i) refers to the methodology used to get various estimates of a given 
quantity, namely i = 1 for the method based on the distribution of first recombination times, 
j = 2 for the cumulative hazard method and i = 3 for the "reactive flux" approach. Tb is the 
average first recombination time and also the typical life time of a chain of average length Lq. 
The estimate r^ is obtained by the long time behaviour of the relevant dynamical function 
used in methodology (i). tv^' is the fraction of transitions which are effective. Such transitions 
can be seen as those which do not belong to sequences of correlated transitions (chain scission 
followed by almost immediate recombination ending into no change in chains topology). An 
additional estimate of tj, (not listed in th etable) may be obtained on the basis of Us and k''' 
computed with methodology i using Th = - — tn — 



LO 


^(«) 


a 


Us * 10^ 


rl^' 


rl^' 


rf^ 


0.1 


.89 


0.89 (0.80) 


0.00759 


694. 


763. 


625. 


0.5 


0.60 


0.68 (0.60) 


0.038 


196. 


212. 


185. 


1.0 


0.45 


0.55 (0.48) 


0.076 


153. 


146. 


122. 


5.0 


0.20 


0.28 (0.25) 


0.38 


60.2 


67.0 


63.3 



Table 1.2. Kinetic data for the semi-dilute case at four scission-recombination attempt fre- 
quencies.All quantities are defined in the previous table 



LO 


^(2) 


^(3) 


a 


n^ * 10^ 


-r* 


rl^' 


-f' 


0.1 


0.85 


0.85 


0.85 (0.73) 


0.0042 


919. 


885. 


813. 


0.5 


0.52 


0.44 


0.53 (0.44) 


0.022 


285. 


265. 


282. 


1.0 


0.34 


0.35 


0.40(0.31) 


0.042 


191. 


190. 


184 


5.0 


0.12 


0.12 


0.16(0.11) 


0.21 


130. 


108. 


106. 



parentheses is the ratio obtained by using data for which the recombination time is 
less than T^ax, while the first value is obtained by assuming that all recombinations 
taking place beyond that maximum time should be necessarily with new partners, 
as it is justified by data shown in figure 1.11 to be discussed later The product aus 
is the number of scissions per unit of time and per unit volume which lead to a 
recombination with a new partner, that is those which lead to real changes in the 
distribution of chain lengths. 

In the explored w range, the value of a indicates that a large fraction of elemen- 
tary scissions are followed by recombinations of the same original chain fragments. 
Using a as a measure of the percentage of effective reactions, we observe a decrease 
of this percentage as the rate of attempted transitions uj increases. Actually, uneffec- 
tive transitions may also come from more complex particular transition sequences. 
Consider a chain end monomer (say monomer i) lying close in space to another chain 
at the level of two adjacent monomers j and k. At high uo, a scission of the bond jk 
may be immediately followed by a recombination between j or k with monomer 
i. In turn, that ik or ij bond may reopen and recombine to restaure the very first 
situation, ending with no effective transition without being detected through the cri- 
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terium of a successive recombination with the same partner. The occurence of such 
events is proven indirectly in the following by noting almost immediate chain end 
recombination with another partner (see figure 1.11). 

To extract estimates of the rate constants defined by the kinetic model of Gates, 
it is thus important to eliminate the spurious transitions from the effective ones. The 
best route for this is certainly to approach the problem via a time scale separation 
between scission-recombination processes taking place on a fast time scale and ef- 
fective transitions taking place on the reaction time scale r^. 

1.6.1 Distribution of first recombination times 

The first approach is to compute the histogram of first recombination times t = 1-2 — ti 
for an arm which became free by scission at time f i and which recombined with an- 
other free arm at time ^2- Figure 1.10 shows the result for the semi-dilute case for the 
four values of tj. In figure 1.11 relative to the case w = 0.5 specifically, we show the 
contribution to the distribution of first recombination times from events implying ei- 
ther the same or a different partner. In the latter figure, one observes that at short time, 
the recombination with the same partner dominates. However, even at short times, 
recombinations with another partner are possible as discussed earlier. The two curves 
corresponding to the two types of contributions cross each other around t = 25 and 
at later times, the recombination with another partner progressively dominates. We 
note however that at times as long as t = 200 where monomers have diffused by a 
distance corresponding to three times their size (see figure 9), there is still 10% of 
recombinations at that time which take place with the same partner. 

Data of figure 1.10 have been analyzed by considering that ^{t) is exponential at 
long times. The slope is interpreted as the inverse of the average life time of a chain 
end which is denoted as r^ . Specific values for the different tj's are indicated in 
tables 1 and 2, respectively for dilute and semi-dilute cases. 

The short time behaviour of the distribution can be analyzed by plotting the func- 
tion versus time on a logarithmic scale. This is done in figure 1.11 where we distin- 
guish contributions from recombinations with the same partner or with a new partner. 
The theory of diffusion controlled recombination kinetics [12] predicts an algebraic 
decay At~^/^ for ^{t) and we observe that it is perfectly satisfied by our data for 
the self recombination part for times larger than t = 1. The picture shows how 
self-recombinations dominate at short times while recombinations with new partners 
dominate beyond t = 30. At very short times (tj 1), a significative part (about 10%) 
of the ^{t) function implies recombinations with anorther partner. 

1.6.2 Cumulative hazard analysis 

To analyse the rates of conformational transitions in butane and other short alkane 
molecules, Helfand [24] suggested to exploit hazard rate plotting. We have adapted 
this technique to the present case. We first summarize the theoretical foundations 
of the method, using the particular case of equilibrium polymer kinetics to directly 
illustrate the concepts. 
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Fig. 1.10. Distribution of first recombination times !f'(t) for the semi-dilute case. Data are 
shown for uj — 0.1 (filled circles), uj — 0.5 (+ symbols), uj — 1.0 (filled lozenges) and 
(jj = 5.0 (filled triangles). Estimates of r^ ' are obtained from the slope of a linear fit of 
In !^(t) versus time in the long time region. 



Let h{t)dt be the probability that a free arm, created at time t = and which 
is still free at time t, undergoes a first recombination in the interval [t,t + dt]. Let 
P{t) be the probability that a free arm, created at time t = 0, has undergone a (first) 
recombination between and t. 

Using these definitions, the following steps can be written 



[1 - P{t + dt)] = [1 - P{t)] (1 - h{t)dt) 
±[l-Pit)] = -hit)[l~P{t)] 



(1.58) 
(1.59) 



[1 - P{t)] = cxp (- f h{t')dt') = exp [-H{t)] (1.60) 

Jo 

where H{t) = L h{t')dt' is the cumulative hazard. In the present case, we anticipate 
a complex process involving correlated events at short times h{i) = h'^^{t) and a 
Poisson process emerging at long times with uniform frequency lirnt^ooh{t) = A 
where A is the recombination rate constant. On the relevant time scale of the kinetics, 
one thus should find for the cumulative hazard function and the P(t) probability 
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Fig. 1.11. Distribution of first recombination times !^(f ) for the semi-dilute case for attempt 
frequencies lo — 0.5 (open triangles). The partition between contributions from recombina- 
tions with the same partner (filled circles) or with a new partner (filled square) are distiguished. 
The power law 0.012(i)^^'"* shows that the "diffusion-controlled" theoretical prediction of 
—5/4 for the power law exponent is perfectly satisfied for the self recombination part for 
times larger than t — 1. The picture shows how self-recombinations dominate at short times 
while recombinations with new partners dominate beyond t = 30. 



H{t) = Hi + Xt 

P{t) = 1 - exp {-Hi) cxp (-Ai) 



(1.61) 
(1.62) 



where Hi, the ordinate intercept of the function H{t) versus t, can be seen as the 
time integral of {h^^{t) — A) from to oo. If a good time scale separation exists, 
eq.(1.62) implies that c = -P(0+) = 1 — exp(— 77/) is the fraction of correlated 
transitions, k'^^ = 1 — c being the estimate according to the present analysis of the 
average probability for a newly created free arm to recombine by the Poisson "mean 
field" kinetic process postulated by Gates in his theory [9]. 

We now explain how the estimate of the cumulative hazard function is con- 
structed. We start by extracting from the BD trajectory a collection of times {£} 
where each member corresponds to the elapsed time between a scission of a particu- 
lar arm at time ti and its next recombination at time t2, so t = t2 — ii. All the 2M 
arms of the system contribute to the t data sample which, for each arm, can contain 
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several times of this kind between the start (at t = 0) and the end (at t — Tmax) of 
the BD trajectory. 

Moreover, if the first change of status of a particular arm since the beginning of 
the BD simulation is a recombination taking place at time t, we can say that this time 
is a lower bound T of an additional unknown elapsed time t between a scission (out 
of our reach) and the next recombination (we observed). Also, if the last change of 
arm status before the end of the BD trajectory is a scission taking place at time t, then, 
the time T = Tmax — t is a similar lower bound of yet another time of interest. The 
analysis thus furnishes a set of K times {t} and I\I lower bounds times {T}) which 
are then separately ordered from the shortest time up to the longer one and indexed 
accordingly as (ti, ^2, ^3, ■•i/<') and (Ti, 12, ^s, ■■■Tm)- Let us consider successively 
all individual arm life times {i}. For any time ti of that collection, the probability 
that a free arm which had survived up to the previous time ii_i changes its status and 
becomes engaged in the formation of new bond between times i^^i and ti is given 
with our available statistics by l/N{ti-i) where the denominator is the number of 
cases (including both types of times t and T) where recombination takes place at a 
time longer than ti_i, i.e. 

7V(tVi) = [K-{i- 1)] + [M - m(,_i)] (1.63) 

where m(j_i) is the index of largest T value which is still inferior to ii_i. In 
terms of these definitions, the cumulative hazard function H{t) can be evaluated at 
each time t 

Hft)^!^ (1.64) 

If the function is linear in time, it implies a Poisson process with a rate of transi- 
tions given by the slope of that linear portion. 

Figure 1.12 shows the resulting cumulative hazard functions for the four lo values 
in the semi-dilute case. A linear behaviour starts around 100 — 150 in all cases and 
it lasts up to times of the order of 1000 where statistical noise sets in. The inverse 
slopes of the linear portions provide estimates t^ ' of the mean life time of a chain 
end. Specific values for both thermodynamic states are gathered in table 1 and table 
2. They are compatible with the values extracted from the histogram of first recom- 
bination times. Using the procedure described at the beginning of this section, one 
gets from the ordinate intercepts of H{t) the n^"^^ estimate indicated in table 2 for 
the semi dilute case. 

The cumulative hazard plot has the advantage that it avoids the necessity of bin- 
ing recombination time data which become scarse at long times. It also exploits sta- 
tistically some additional information (useful for low rates) from portions of the 2M 
trajectories between a scission and the next recombination in a way which is trun- 
cated either at the beginning or at the end of the BD run. 
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Fig. 1.12. Cumulative Hazard curves of recombination times for the semi-dilute case. Data are 
shown for various attempt frequencies, namely co — 0.1 (dashed line), co — 0.5 (dot-dashed 
Une), uj = 1.0 (dotted line) and u = 5.0 (continuous line). 



1.6.3 The reactive flux correlation function and the transmission coefficient 



We now discuss still another way to isolate short time transitional effects in order to 
estimate the rate of "effective" transtions. It relies on the way so-called "Transition 
State Theory" (TST) values of chemical reactions rates can be corrected by a mul- 
tiplicative transmission coefficient < k < 1 which takes into account the fraction 
of events where the products formed by the forward reaction are effective transi- 
tions. By effective, it means that after such a transition, the product is able to relax 
its energy excess (which was accumulated to cross the forward barrier) sufficiently 
quickly to avoid an immediate recrossing of the barrier in the backward direction. 

In our case, the analogy is clear. If a bond scission is followed by the same bond 
recombination within a short microscopic time scale (too short for relative diffusion 
to take place), it basically means that no effective scission has happened to the chain 
to which the bond belongs on time scales larger than the relaxation time of a chain 
fragment of a few monomers. 

We have thus estimated k'-^^ by exploiting the concept of the reactive flux correla- 
tion function, as explained by Chandler in the framework of an equilibrium between 
forward and backward unimolecular reactions [19]. We assume here that for an equi- 
librium between forward unimolecular reactions (chain breaking) and bimolecular 
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backward reactions (end-chain recombinations), the same conceptual expression can 
still be applied for the forward reaction. 

Let us consider one specific arm (among the two of one particular monomer) 
which is labelled by an index i and to each arm, we associate a signature variable 
Si {t) which is equal to unity if the ith arm is free at time t and equal to zero if the 
arm is engaged in a binding with another monomer at that time. Crossings are then 
characterized by a step change of the Si{t) values, Si{t) going from to 1 or from 
1 to for scissions and recombinations respectively. A statistics over all crossing 
events is defined: for each individual crossing, we set the time at crossing to zero, 
denoting times very shortly after or before the transition by t = 0+ and t = 0~. The 
reactive flux correlation function is then defined by[25] 

C'RFit) = 2 < [5,(0+) - 3,(0-)] S,it) >,ross (1.65) 

where the average is performed over all crossings (both scissions and recombi- 
nations). If we treat separately the averages on scissions and recombinations, we get 
the equivalent alternative formulation 

^RFy^) ^ 'Jiy^j '^scissions ^ ^i\^) ^recombinations (,i.uOj 

This function starts from unity at very short (positive) times with the first term 
being one and the second being zero. On very long times (long respect to the relevant 
inverse kinetic rate constants), this function goes to zero as it is the difference of two 
averages which ultimately converge to the same equilibrium average value < Si >. 

If all transitions are effective, the reactive flux function will decrease smoothly 
from 1 to zero reflecting only the kinetics of these effective recombination and scis- 
sions. If transitions are often followed by opposite transitions within a short time 
scale Tynicro or if successive transitions are correlated within such a short time, this 
will be reflected in the reactive flux by a fast initial decay over that microscopic time. 
At later times, the macrocopic kinetics (the one considerd by mean field theoretical 
expresions) will be reflected by the Crf (t) evolution towards zero on a macroscopic 
time scale, Tmacro- The transmission coefficient k^'^^ is thus estimated as the initial 
value of the "slow" decay part of that Crf (t) function. This can be expressed as 

K^^^ = CRF{T^^cro) (1-67) 

with the requirement of a good time scale separation, that is Tmicro "C Tmacro- 

Figure 1.13 shows the u = 0.5 case for the semi-dilute case already discussed 
in fig 1 . 11 in the context of the distribution of first recombination after a scission. 
This curve is more difficult to analyse because it is difficult on the basis of this sole 
function to decide when the behaviour of ln{CRF{t)) really becomes close from a 
linear behaviour We had to rely on our experience with the distribution 9{t) and the 
cumulative hazard function H{t) for that state point and w value. As both functions 
suggest that the simple "mean field" process emerges at times later than 150, we have 
fitted the CRF{t) data by an exponential function Cj^p{t) = k^^^ * exp {—t/rl ) in 
the accessible time window above t = 150. In the present case, we see that there 
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Fig. 1.13. Reactive flux correlation function in the semi-dilute case for uj = 0.5. The dashed 
line corresponds to a fit by a function «: exp {—t/rt) in the time region 150 < t < 500. Best 



parameters found are k 



(3) 



: 0.48 and t. 



(3) 



282. 



is no real time scale separation between a short time behaviour, observed in the 

(3) 

< < < 150 regime, and the macroscopic time scale r^^ ' which is equal to w 280. 
Nevertheless, the value of r^^' ' turns out to be rather close to the other estimates of a 
chain end average life time and the k^'^' value will be shown in the next subsection 
to lead to a reasonable correction of Ug, the total number of transitions per unit of 
time and unit volume quoted in tables 1 and 2. Results of k^^^ and r^ for all uj's 
and both state points are included in the tables 1 and 2. 



1.6.4 Estimation of rate constants: comparison of the various methods 

Explicitly, the scission rate constants kg and kr of the "mean-field" kinetic model 
can be estimated from our data in two distinct ways. 

• On the basis of the total number of transitions Ug per unit time and per unit 
volume, and on the basis of an estimate of k by one of the methods discussed 
above, the rates can be estimated from a trivial modification of eqs. (1.40,1.41), 
namely 
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KUs = ks4> = -^—^ (1.68) 

• On the basis of the long time behaviour of a dynamical function relevant to the 
particular methodology adopted, the "chain end" life time can be directly esti- 
mated. As the latter is equivalent to the typical life time of a chain of average 
length denoted as t;,, one gets rate constants through the equivalence 

1 2io 

KgLf) krCp 

where the last equality follows from detailed balance requirements. 

Tables 1 and 2 provide all needed data. We find that the implicit time scale sep- 
aration inherent to all analyses provide compatible estimates of the effective rate 
constants. This can be verified by applying eqs (1.68) and (1.69) with the particular 
estimate of the left hand side. 

Looking at all values in the tables, we get an overall consistency between all 
methodologies. It indicates that all our strategies to extract the "macroscopic" rate 
constants do work. Among them, the cumulative hazard analysis appears the most 
straightforward to analyse once H{t) is known and particularly robust as noise influ- 
ence seems minimal. 

We are now in a position to compare the percentage of effective transitions n 
based on a time scale separation and the a value based on the percentage of recom- 
binations involving a new partner respect to the one left by scission. We observe in 
the two tables that both estimates are reasonably consistent and quickly decrease for 
io increasing. Fine analysis suggests however that some distinction exists between 
both quantities in the present case, that is a free Rouse relaxation type of dynamics 
in solution for moderately long chains. It requires further analysis and the consid- 
eration of more entangled cases to make progress on that point. Finally, it should 
be stressed that while diffusion controlled and mean field kinetics have sometimes 
been opposed to each other[16], it appears in our analysis that they consistently de- 
scribe short time correlation effects dominated by self recombinations and long time 
kinetics dominated by recombinations with new partners. 

1.6.5 Analysis of the monomer diffusion 

The mean squared displacement (MSD) of all monomers as a function of time is 
shown in Figure 1.14 for all cj's in the semi-dilute case. As these monomers belong 
to a polydisperse set of chains which continuously break and recombine, the inter- 
pretation of the MSD is rather complex. 

• Being a single monomer property, a relevant quantity is the average fraction of 
monomers belonging to chains of length L, namely a distribution oc L * co{L). 
Fig 1.15 shows the resulting curve for the adopted semi-dilute conditions with 
<L >=49.9. 
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Fig. 1.14. Mean square displacement of monomers (semi-dilute case). The continuous curve 
has been obtained for the relevant polydisperse sample but with u set to zero, averaging over 
many independent initial configurations. It is a reference from which the other curves depart 
at decreasing times when increasing values of uj are considered. The cases lj = 0.1 (dashed 
lines), w = 0.5 (dot-dashed lines) and u — 5.0 (dotted lines) are shown together with two 
indicative power laws which suggest the evolution of the subdiffusive scaling regime to the 
normal diffusion limit. 



The time scale for the scission/recombination process has been estimated so far in 
termsof the average survival time l/(fcs < L >) for the polymer of average size, 
a quantity which was found to vary from w 100 up to w 900 in the investigated 
u! range (see table 2). In the present context, we need to consider the explicit L 
dependent average survival time oc l/{ksL) for a given cu value and figure 1.15 
shows the resulting functions for two such frequencies. 

To discuss the monomer mean squared displacement, it is also important to have 
an estimate of the longest internal relaxation (Rouse) time of the chains. Perform- 
ing a separate brownian dynamics run on a monodisperse sample of 20 "dead" 
polymers of length L=50 (corresponding to the mean polymer length of our liv- 
ing polymer sample) at the same semi-dilute state point (T — 1, (f) ~ 0.15), we 
get from the long time exponential bevaviour of the end-to-end vector relaxation 
function a value of r/j ~ 4.24 10^. For our very flexible polymers, the scaling of 
the Rouse time with L in the semi-dilute conditions should scale as 
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tr{L)^AL\ (1.70) 

On the basis of our estimated value for L = 50 chains, we get Ak, 1.70, a value 
which can be compared (at least its order of magnitude) to theoretical estimate 
[6, 26] A K, Tuob/{L*Y ™ terms of Tuob, the blob relaxation time, and L* , 
the number of monomers per blob. We get A = 5.8 for our conditions (see 
discussion in the section on chain conformational properties) where the number 

3 

of monomers per blob L* « 37 and Tuob ~ ^^^^ ~ 8000 using a blob size 
X = 20 and a solvent viscosity 7]s = -^^ = 1 estimated from the friction 
coefficient using the Stokes law with an hydrodynamic radius fixed to a/2. 

On the basis of the above considerations, it is useful to define a particular 
chain length A for which the Rouse time is equal to the survival time, giving 
A = (1.7 * fc^)-^/^, that is yl = 30, 20, 18, 15 for uj = 0.1, 0.5, 1.0, 5.0 respec- 
tively. Chains longer than A have their internal dynamics strongly altered by scis- 
sions, while chains shorter than A should have an internal dynamics little affected 
with respect to dead polymers of the same size. As figure 1.15 illustrates, when 
the scission attempt frequency increases, the fraction of monomers which belong to 
chains whose dynamics is little affected by scissions decreases. 

Up to a time of 100, the MSD shown in figure 1.14 is uj independent and presents 
a Rouse like behaviour with a time scaling exponent t°^. In the reference case where 
u) was set to zero (we simulate polydisperse sample of dead polymers), this power 
law bahaviour persists at least up to time 10000, a time corresponding to the Rouse 
time of a chain of length 170. When scissions and recombinations are allowed, the 
MDS deviates from the master curve (of dead polymers) to adopt progressively a 
MSD linear in time. This takes place sooner and sooner for increasing values of uj. 
At long times (t > 2000), the MSD has evolved towards usual Einstein diffusion 
with a diffusion coefficient increasing with kg. In figure 1.16, D is observed to (very 
roughly) scale like fc^'^^ in the explored range, a result which is not so far from 
the 1/3 power behaviour observed by Milchev with the lattice BFM lattice [14], 
but which requires further analysis with better estimates of D as statistics is poor 
above t=10000. The 1/3 power law scaling is explained by assuming that clusters of 
A monomers are responsible of the long time behaviour of living polymers, giving 

1.6.6 Analysis of the first recombination times distribution and corresponding 
diffusive steps 

In figure 1.17, we report for two different bond change attempt frequencies cj, the 
mean squared displacement of a chain end between its creation and its next recombi- 
nation, as a function of the corresponding life time of that chain end. Interestingly, we 
observe a common behaviour of this dependence for both frequencies, meaning that 
the diffusive dynamics of chain ends is little sensitive to the scission/recombination 
processes of the chains as a whole. 
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Fig. 1.15. Characteristic times as a function of polymer size. The Rouse relaxation time tr = 
1.70*1/^, shown as dotted curve, must be compared to the lifetime of a polymer of size L, as 
given Tb = {ksL)~^ . Two frequencies are illustrated, namely uj — 1.0 (dot-dashed line) and 
a; = 0.1 (dashed line) using effective ks values. We show on the same graph the distribution 
of monomers (continuous curve) as a function of the size of the polymer they belong to (see 
text). 



As chain ends have very different life times for different tj's, the typical average 
distance travelled by a chain end between its creation and its recombination is also 
function of that frequency. Mean travel distances of 3.8, 4.6 and 6.8 are obtained for 
the average chain end life times at lo =5.0, 1.0 and 0.1 respectively, to be compared 
with the mean static distance between chain ends (assuming no spatial coiTelations) 
given hy h ~ {< L > /20)3 which turns to be /i = 5.5 in our semi-dilute case. 
By comparison, the radius of gyration of the average chain length in the system is 
Rg{< L >) K, IQ. The dynamic picture is thus that chain ends diffuse similarly with 
time in an environment of other chains which continuously break and recombine. As 
the reacting frequency increases, the probability of a chain end to recombine with a 
neighbouring monomer increases and thus takes place earlier, after diffusion over a 
shorter distance. 
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Fig. 1.16. Monomeric Einstein diffusion coefficient as a function of ks for the semi-dilute 
case. The value of D was obtained from the long time behaviour (beyond t=8000) of the 
monomeric mean square displacement shown in figure 1.14. The dashed line represent a k^'^^ 
power law. 



1.7 Summary and conclusions 



In this paper, we have proposed a new version of a mesoscopic model of equilibrium 
polymers in solution and we have exploited this model to study the structure and the 
dynamics of these complex systems at a dilute and at a semi-dilute state point. 

The dynamics at the semi-dilute state point has been more specifically investi- 
gated. Our system consists of chains with an average number of Lq « 50 monomers 
which can be viewed as a collection of blobs of « 37 monomers. Given the large 
flexibility of the model and the moderately semi-dilute character, our polydisperse 
polymer solution should be in a non-entangled dynamic regime, with a chain relax- 
ation time growing with chain length L as L^. 

Working with a unique chain length distribution governed by the structural and 
thermodynamic microscopic parameters, we have analyzed the specific influence 
of the Monte Carlo binding/unbinding change attempt frequency on the dynami- 
cal properties. Mean field theories are found to be valid for times of the order and 
beyond the mean life time of a chain of average size provided the kinetic constants 
computed from the Monte Carlo accepted binding/unbinding changes are rescaled by 
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Fig. 1.17. Mean square displacement of end chain monomers between a time set to zero when 
produced by scission up to a time t corresponding to their first recombination, as a function 
of their total life time t. Due to poor statistics, only data for to — 1 (filled squares) and uj — 5 
(empty circles) are shown and are limited to t = 10000 maximum. Data on the mean square 
displacement of all monomers for corresponding uj values already shown in figure 1.14 are 
shown again here for comparison (see text) up to a maximum time of 20000. Scaling power 
laws l.Ot^" (thick dashed line) and 0.943t°'^ (thick dot-dashed line) shown up to t=20000 
are used in the text. 



a transmission coefficient which is a measure of the fraction of successful scissions. 
This transmission coefficient, when estimated on the basis of a time scale separa- 
tion, turns out to be close to the fraction of recombinations involving the reunion of 
a chain end with a new partner and no longer with the monomer it was attached to 
previously. 

The rate constants have been estimated by different techniques which give co- 
herent estimates. The cumulative hazard function introduced by Helfand to compute 
isomerisation rates in chain molecules has been found to be particularly efficient and 
simple to implement. In paticular, we find agreement between 1) kinetic constant 
estimates based on the the total number of transitions rescaled by the transmission 
coefficient coming from a diffusion controlled mechanism and 2) estimates based 
directly on the time scale of the long time exponential decay of relaxation functions 
(mean field Poisson process). It would be interesting to perform so-called T-jump 
experiments [26] with our system to test whether the time evolution of the average 
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polymer length Lq towards its new equilibrium value at a different temperature con- 
firms the mean field analytical prediction with the kinetic constants computed in the 
present work from various population relaxation functions at equilibrium. 

The chain size that we chose originally was adequate to demonstrate the validity 
of the methodological aspects developped in this work. In order to test dynamical 
scaling with a wider range of the scission-recombination process frequency, longer 
chains will be needed which may require better optimized methods to take solvent 
effects into account. The extension of this study towards the rheological properties 
of our model system are under way. 
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